home *** CD-ROM | disk | FTP | other *** search
/ MacFormat 1994 August / August CD.bin / Shareware / Education / numericalmethods Folder / chap_11 / jacobi2.m < prev    next >
Encoding:
Text File  |  1994-06-05  |  1.3 KB  |  46 lines  |  [MATF/MATL]
  1. function [V,D] = jacobi2(A,epsilon,show)
  2. % [V,D] = jacobi2(A,epsilon,show)
  3. % To compute the eigenpairs of a symmetric matrix.
  4. % The cyclic Jacobi`s method of iteration is employed.
  5. % A is an n x n matrix, input.
  6. % epsilon is the tolerance, input.
  7. % V is the diagonal matrix of eigenvectors, output.
  8. % D is the diagonal matrix of eigenvalues,  output.
  9. if nargin==2, show = 0; end
  10. D = A;
  11. [n,n] = size(A);
  12. V = eye(n);
  13. cnts = 0;
  14. cntr = 0;
  15. done = 0;
  16. working = 1;
  17. stat = working;
  18. while (stat==working),
  19.      cnts = cnts+1;
  20.      t = sum(diag(D));
  21.      stat = done;
  22.     for p = 1:(n-1),
  23.       for q = (p+1):n,
  24.            if ((abs(D(p,q))/t)>epsilon),
  25.           t = D(p,q)/(D(q,q) - D(p,p));
  26.           c = 1/sqrt(t*t+1);
  27.           s = c*t;
  28.           R = [c s; -s c];
  29.           D([p q],:) = R'*D([p q],:);
  30.           D(:,[p q]) = D(:,[p q])*R;
  31.           V(:,[p q]) = V(:,[p q])*R;
  32.                 cntr = cntr+1;
  33.           if show==1,
  34.             home; if cntr==1, clc; end; 
  35.             disp(['Jacobi iteration No. ',int2str(cntr)]),disp(''),...
  36.             disp(['Zeroed out the element  D(',num2str(p),...
  37.                   ',',num2str(q),') = ']),disp(D(p,q)),...
  38.                disp('New transformed matrix  D = '),disp(D)
  39.           end
  40.                 stat = working;
  41.         end
  42.       end
  43.     end
  44. end
  45. D = diag(diag(D));
  46.